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A  Ihree^menslonal  ItMoiy  Is  derived  to  describe  the  temporal  behavior  of 
gravity  curtenls  of  cohestonless  gtorwlor  media  In  on  attantof  to  model  the 
motion  of  dense,  flow-type  snow  avalanches.  Ice  and  rock  slides.  A  Mohr- 
Ooutomb  yIeU  criterion  Is  assurned  to  describe  the  constitutive  behavior  of  the 
material,  and  Ihe  basal  bed  fttcHon  Is  described  sknltorly  by  0  Coulomb  type  of 
frtctton.  Adrog  term  Is  Included  In  ordertotnodel  the  occurrence  of  flow  regltnes 
where  boundary  drag  becomes  non-negllgible.  Data  from  laboratory  slmulo- 
ttons  ore  compared  to  a  series  of  numericol  studies  based  on  the  aforemen¬ 
tioned  theory.  A  rrandlmenstonol,  depth  and  width  averaged  form  of  the  theory 
Is  considered.  A  Logronglan  Unite  difference  scheme  Is  then  applied  to  numeri¬ 
cally  model  some  limiting  coses  of  the  governing  equations.  Two  different 
numerical  models  are  developed,  tested  and  compared  to  experimental  values. 
The  results  Indicate  that  the  model  can  account  for  flow  transitions  by  Inclusion 
ofttw  drag  term  vt^ien  the  Initial  Inclination  angle  Is  large  enough  to  affect  bound¬ 
ary  drag.  Furthermore,  the  temporal  and  spatial  evolution  of  the  granulate  and 
final  runout  position  con  be  predicted  to  values  well  within  the  experimental 
error. 
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NOMENCLATURE 


2*  Cartesian  coordinates 

z*  curvilinear  coordinates 

natural  basis  for  curvilinear  system 
Qi  j  metric  of  the  space 

orthonormal  basis  for  curvilinear  system 

V  underline  denotes  first  order  tensor  (vector) 

t  double  underline  denotes  second  order  tensor 

C  angle  defining  longitudinal  variations  from  the  horizontal 

u<i>  physical  component  of  v 

p  mass  density 

V  del  operator,  operating  on  the  argument  indicated  by  arrow 

f <i  j>  physical  component  of  for  the  stress  tensor  with  the  subscripts  implying  on  the  i  surface 

and  in  the  j  direction 
/  function  defining  free  surface 

b  function  defining  bed  surface 

n  unit  normal  vector 

V  velocity  vector 

k  body  force  vector 

a  acceleration  vector 

t  Cauchy  stress  tensor 

6  bed  friction  angle 

(j>  internal  friction  angle 

Vy  boundary  drag  coefficient 

£*  characteristic  length  scale 

Q  acceleration  due  to  gravity 

T  characteristic  time  scale 

X  overbar  denotes  nondimensional  quantity 

X  double  overbar  denotes  nondimensional,  depth-averaged  quantity 

X  tilde  denotes  nondimensional,  depth-  and  width-averaged  quantity 

6t  i  ratio  of  £’  characteristic  length  scale)  to  P  characteristic 

length  scale) 

tx7  Voellmy  coefficient 
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Model  for  Avalanches  in  Three  Spatial  Dimensions 

Comparison  of  Theory  to  Experiments 

RENEE  M.  LANG  AND  BRIAN  R.  LEO 


INTRODUCTION 

Increases  in  land  use  and  development  in  mountainous  areas,  anticipated  warming  of  the 
earth’s  atmosphere,  and  the  deforestation  in  alpine  regions  from  air  pollution,  acid  rain, 
previously  poor  forestry  practices,  and  land  misuse  have  renewed  attention  to  the  impor¬ 
tance  of  reliably  recognizing  natural  avalanche  paths  and  predicting  deposition  zones.  The 
impact  pressures  induced  by  an  avalanche  in  the  runout  area  motivate  zoning  restrictions 
and  building  code  requirements.  Existing  avalanche  models  are  generally  inadequate.  There¬ 
fore,  the  current  trend  in  the  research  of  flowing  snow  and  other  materials  is  still  directed  at 
determining  the  predominant  mechanisms  governing  the  motion. 

The  phenomena  collectively  referred  to  as  avalanches  can  be  physically  characterized  as 
multiphase  gravity  flows,  which  consist  of  randomly  dispersed,  interacting  phases,  whose 
properties  change  with  respect  to  both  time  and  space.  Lang  and  Dent  (1982)  describe  an 
avalanche  as  “the  transient,  three-dimensional  motion  of  a  variable  meiss  system  made  up 
of  nonrigid,  nonrotund,  nonuniform  assemblage  of  granular  (snow)  fragments  flowing  down 
a  nonuniform  slope  of  varying  surface  resistance.”  In  this  sense,  an  exact  analysis  of  the 
motion  of  an  avalanche  is  perhaps  an  unattainable  goal. 

Currently,  the  most  widely  used  and  applied  avalanche  models  utilize  a  center  of  mass 
approach  (i.e.,  the  mass  is  accretized  at  one  point)  and  are  based  on  ideas  first  suggested 
by  Voellmy  (1955).  By  using  the  Prandtl  mixing  length  approach  (see  for  example  Pao 
1961),  Voellmy  related  the  stress  at  the  base  of  the  flow  to  the  square  of  the  velocity,  and 
then  postulated  an  additional  Coulomb  friction  force.  However,  he  assumed  uniform,  steady 
conditions  and  derived  a  method  that  is  diflicult  to  apply,  in  the  sense  that  a  number  of 
subjective  parameters  must  be  predetermined  in  order  to  obtain  results  that  match  observa¬ 
tional  data.  Although  Voellmy’s  model  has  been  successively  improved  (Salm  1966,  Mellor 
1968,  Perla  et  al.  1980),  these  models  have  not  advanced  beyond  the  aforementioned  center 
of  mass  approach.  The  height  of  the  flow  may  be  included  as  a  parametric  value,  but  is  not 
calculated  as  a  function  of  space  and  time. 


Other  models  attempt  to  idealize  these  complicated  materials  as  linear  Newtonian  fluids 
(Brugnot  1979,  Dent  and  Lang  1980),  leading  to  the  Navier-Stokes  equations  that  may  be 
solved  numerically.  Although  this  type  of  constitutive  relation  perhaps  should  not  be  assumed 
to  adequately  describes  these  media,  some  success  has  been  achieved  in  modeling  certain 
aspects  of  avalanche  behavior  with  this  approach.  Apparent  viscosities  may  be  measured  for 
fluidized  solids  (e.g.  Siemes  and  Hellmer,  1962).  The  ease  with  that  particles  fluidize  and  the 
range  of  conditions  that  sustain  fluidization  vary  greatly  among  solid-fluid  gravity  current 
systems. 

When  considering  snow  avalanches  in  particular,  due  to  the  nonlinear  behavior,  and 
strong  stress  and  strain  history  dependence  of  snow,  constitutive  equations  tend  to  be  quite 
complicated.  Both  continuum  theories  (Salm  1971,  1975,  1977,  Brown  et  al.  1973,  Brown 
and  Lang  1973,  Brown  1976,  1977,  1979,  1980a,  1980b)  and  microstructural  theories  (St. 
Lawrence  1977,  St.  Lawrence  and  Lang  1981,  Hansen  and  Brown  1986,  Hansen  1985)  have 
been  proposed  and  applied  successfully  to  various,  specific  problems  in  snow  mechanics.  But 
to  some  extent  all  require  detailed  information  on  the  particular  snow  type  in  question,  that 
entails  a  wide  range  of  possible  parametric  values. 

Another  hindering  factor  has  been  the  availability  of  experimental  data  on  avalanche 
motion.  Field  measurements  are  limited  and  difficult  to  obtain  (see,  for  example,  Norem  et 
al.  1985,  Gubler  1987),  and  similitude  is  difficult  to  achieve  when  attempting  to  scale-down 
and  model  such  flows  in  the  laboratory.  For  example,  Lang  and  Dent  (1980)  utilized  a  1/100 
scale  model  of  snow,  and  successfully  matched  the  Froude  number,  but  not  the  Reynolds 
number.  This  is  not  an  unusual  problem  in  scale  modelling  when  the  same  “fluid”  is  used. 
For  a  scale  model  with  length  1/n  of  the  prototype,  equality  of  the  Reynolds  number  requires 
a  model  characteristic  speed  of  Vm  =  nvp,  where  Vp  is  the  prototype  characteristic  speed. 
Equality  of  the  Froude  number  requires  Vm  =  Vp/  y/n.  Both  requirements  cannot  be  satisfied 
by  using  the  same  “fluid”  in  both  model  and  prototype.  Other  recent  attempts  to  scale  model 
avalanche  flow  include  Nakamura  et  al.  (1987),  Pliiss  (1987),  Hutter  et  al.  (1988),  and  Lang 
et  al.  (1989). 

Savage  and  Hutter  (1989)  demonstrated  that  one  may  adequately  describe  a  flowing 
mass  of  cohesionless,  granular  medium  as  an  incompressible  continuum  that  obeys  a  Mohr- 
Coulomb  yield  criterion.  If  the  moving  mass  is  assumed  to  be  long  and  thin,  a  nondi- 
mensionalization  scheme  (resembling  shallow-water  theory)  may  be  employed,  which  renders 
more  tractable  governing  equations.  By  neglecting  the  boundary  layer,  a  classical  approx¬ 
imation  method  (first  introduced  by  von  Karman  (1921)  and  later  refined  by  Pohlhausen 
(1921))  is  employed,  in  that  the  equations  are  averaged  over  depth  and  a  velocity  profile  is 
assumed.  By  defining  a  depth-averaged  streamwise  velocity  and  eissuming  incompressibility, 
these  equations  may  be  reliably  solved  with  a  Lagrangian  finite  difference  scheme,  such  that 
the  computational  grid  is  advected  with  the  material. 

The  purpose  of  this  investigation  was  to  extend  this  theory  to  three  spatial  dimensions  and 
allow  for  the  inclusion  of  a  boundary  drag  term,  and  test  its  applicability  against  experimental 
results.  The  governing  equations  are  derived  in  a  curvilinear  coordinate  system  that  allows 
for  topography  to  vary  as  in  a  natural  avalanche  path.  The  only  assumptions  included  in  the 
development  are  incompressibility  and  the  application  of  a  Mohr-Coulomb  type  material  law, 
while  also  allowing  the  inclusion  of  a  boundary  drag  term.  A  nondimensionalization  scheme 
is  employed,  choosing  three  characteristic  length  scales.  By  depth  and  width  averaging  the 
balance  equations  for  the  theory  in  three  spatial  dimensions,  sidewise  spreading  of  the  mass 
may  be  included.  A  Lagrangian  finite  difference  scheme  is  applied  to  numerically  integrate 
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the  sets  ot  equations  that  are  derived.  The  numerical  simulations  are  then  compared  to  the 
laboratory  findings.  The  experimental  results  have  been  presented  in  Lang  et  al.  (1989). 

The  analytical  results  indicate  that  a  flow  transition  regime  exists  that  is  analogous  to  a 
laminar  to  turbulent  transition  in  a  fluid,  and  determines  the  applicability  of  the  boundary 
drag  term.  Results  are  reviewed  and  suggestions  are  made  for  appropriate  investigations 
implied  by  these  results. 


MATHEMATICAL  FORMULATION 

The  governing  equations  and  closure  relations  are  proposed  to  represent  the  flow  of  a 
cohesionless  granular  medium  on  any  given  topography.  By  describing  granular  media  as 
multiphase  materials,  intraphase  balance  equations  for  mass,  momentum,  and  energy  and  the 
associated  boundary  and  initial  conditions  may.  in  principle,  be  written.  But  the  resulting  set 
of  equations  is  usually  too  complicated  to  allow  for  a  solution.  For  example,  the  application  of 
kinetic  theories  for  granular  flow  would  in  general  require  the  solution  of  the  energy  equation, 
determination  of  granular  temperatures,  and  velocity  and  density  gradients.  Szidarovsky  et 
al.  (1987)  discuss  the  inherent  difficulties  in  attempting  to  apply  these  theories  to  chute  flows 
of  granulate.  Furthermore,  such  information  is  inappropriate  for  an  engineering  application. 

The  goal  is  therefore  to  develop  a  simple  model  to  predict  the  flow  characteristics  of  the 
mass  of  granulate.  A  purely  mechanical  theory  is  proposed  and  the  granular  mass  is  treated 
as  a  continuum.  The  concept  of  representing  the  granulate  as  a  continuum  should  pose 
no  conceptual  problem,  as  a  mixture  theory  approach  would  require  the  sum  of  the  phase 
interaction  terms  to  vanish  in  the  overall  balance  equations,  in  any  case. 

The  requisite  is  a  suitable  constitutive  assumption.  In  order  to  justify  the  following 
approach,  it  is  necessary  to  consider  the  different  flow  regimes  in  question,  and  also  the 
experimental  results  from  other  investigations. 

First,  two  extreme  flow  regimes  for  granular  media  are  commonly  recognized  (see,  for 
example,  Savage  1984).  At  very  low  shear  rates,  a  granulate  may  be  characterized  as  a 
Mohr-Coulomb  material.  This  implies  quasistatic  deformation,  to  some  degree.  The  Mohr- 
Coulomb  law  is  the  basis  of  Critical  State  Soil  Mechanics  (CSSM),  and  is  commonly  referred 
to  as  the  rupture  condition;  i.e.,  it  implies  impending  plastic  flow  for  a  static  soil  mass,  as 
it  predicts  the  limit  equilibrium  when  the  shear  stress  is  equal  to  the  shear  strength  of  the 
noncohesive  soil.  Simply  stated,  the  maximum  ratio  of  the  shear  to  normal  stress  is  equated 
to  the  tangent  of  the  internal  friction  angle  of  the  granulate.  For  slowly  flowing  masses  of 
granulate,  the  bed  shear  stress  results  from  the  dry  friction  forces  at  the  particle  contacts, 
and  the  Coulomb  failure  criterion  adequately  describes  these  types  of  flows. 

When  shear  rates  are  very  high  and  particles  are  widely  dispersed,  instantaneous  collisions 
dominate  the  behavior  of  the  granulate.  Interparticle  and  particle-bed  contacts  are  brief, 
and  the  fluctuating  part  of  the  particle  motion  cannot  be  neglected.  Kinetic  theories  predict 
quadratic  dependence  of  stress  on  shear-rate,  for  rapidly  shearing  highly  dispersed  particulate, 
and  this  seems  to  provide  good  agreement  with  the  available  laboratory  data.  This  extreme 
case  may  be  used  to  describe  the  behavior  of  the  dust  cloud,  which  only  serves  to  obscure 
the  true  nature  of  the  avalanche,  but  does  not  represent  any  type  of  snow  avalanche  activity 
that  could  be  construed  as  destructive,  and  need  not  be  considered. 

The  real  avalanche  flow  situation  in  most  cases  may  develop  from  a  low  shear  rate  flow  into 
an  intermediate  regime  where  both  long-term  particle  contact  and  instantaneous  interparticle 
collisions  contribute  to  the  generation  of  stresses.  This  was  shown  to  be  the  case  in  the 
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previously  discussed  experiments  in  Lang  et  al.  (1989),  for  the  artificial  avalanches  released  at 
angles  exceeding  35".  Observations  (Dent  and  Lang  1981)  and  Doppler  radar  mesisurements 
(Gubler  1987)  on  natural  snow  avalanches  confirm  that  a  dense  slowly  deforming  mass  rides 
on  a  rapidly  shearing  layer,  and  may  or  may  not  be  encompassed  by  a  highly  suspended  dust 
cloud.  This  type  of  flow  is  more  difficult  to  handle  theoretically,  as  the  exact  mechanisms  have 
not  yet  been  clearly  defined.  Savage  (1983)  suggested  representing  the  total  stresses  as  a  linear 
sum  of  the  rate-independent  dry  friction  and  a  rate-dependent  “viscous”  part.  The  question 
remains  as  to  how  the  dynamic  friction  angle  actually  changes  at  lower  volume  fractions  and 
higher  shear  rates.  Following  the  work  of  Savage  and  Hutter  (1989),  the  governing  equations 
of  mass  and  momentum  balance  are  derived  in  three  spatial  dimensions  in  a  coordinate  system 
that  allows  for  any  given  topography.  An  angle  (  is  defined  such  that  C  is  a  function  of  x*, 
the  longitudinal  direction,  and  defines  the  local  inclination  angle  of  the  generatrix, 

C  =  C(^*)  •  (1) 

It  is  then  possible  to  define  Cartesian  components in  terms  of  the  curvilinear  coordinates  x* 
and  the  natural  basis  {  5^  }■  An  orthonormal  basis  (  ^  }  in  terms  of  the  covariant  components 
of  the  metric  or  fundamental  tensor  of  the  space  Qmn  =  £„,  •  £„  is  defined  and  applied. 

The  Eulerian  form  for  mass  balance  is 

-t-p(x,t)  •  t;(x,i)  =  0  .  (2) 

In  terms  of  the  physical  components  of  v  and  assuming  that  the  granulate-air  mixture  is 
incompressible,  mass  balance  is  written  as 


dv<\>  dv<2>  dv<z> 

»  r\  o  ' 


dx^ 


dx"^ 


dx^ 


I  ^.^3  ( ,  *!<3>  ^  _  n 

5x1  [  Qx2  +  53.3  dx^  ® 

The  principle  of  balance  of  linear  momentum  can  be  written  as 


(3) 


^•t  +  pk  =  P&  ^  (4) 

where  t  is  the  Cauchy  stress  tensor,  6  is  the  body  force  vector  and  a  is  the  acceleration 
vector.  The  equations  of  motion  can  be  expressed  in  terms  of  the  physical  components  as 
follows: 

in  the  xi  length  direction  as 


^^<1 1> 
5x1 


=  P 


9x' 
5v<i> 


5<<2  1>  ,  5t<3  1>  ^ 

+  -t-pgsmC 


5x3 


5x3 


dv<i> 
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dv<i>  ,  5v<i> 
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5x3 


(5) 


in  the  x^  width  direction  as 
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(6) 
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\  fdt<2 2> 
)  .  dx'^ 


+ 


9t<3 2> 
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dx^ 
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dt 
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+  V<3> 
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and  in  the  x^  depth  direction  as 


9t<l3> 

dx^ 


+  ^(^<3  3>  -  tci  1> 


dt<2 
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dx"^ 
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dt 
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dx^ 
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\  at)<3> 
/  dt 


+  u<2> 


^<3> 

dx^ 


+  W<3> 


^^<3> 

dx^ 


(7) 


The  above  forms  for  the  balance  of  momentum  equations  are  shown  in  order  to  emphasize 
the  coupling  of  many  terms  with  the  degree  ci  curvature  dC,/dx^  and  the  depth  of  the  flow 
x^. 

The  boundary  conditions  need  to  be  defined  at  the  solid  surface,  or  base  of  the  flow,  and 
at  the  fi:ee  surface. 

By  defining  the  free  surface  as 


f  =  f{x^,x^,t)  =  x^  ,  (8) 

and  requiring  the  local  and  convected  rates  of  change  of  /  to  balance,  the  well-known  kine¬ 
matic  condition  is  obtained. 


■^+V/  v  =  0  ,  (9) 

which  simply  states  that  the  firee  surface  is  material,  and  applies  when  /  is  generated  by 
streamlines. 

In  terms  of  the  physical  components,  the  kinematic  condition  reads 


\df  ^  df 


=  0 


(10) 


A  balance  of  forces  at  the  free  surface  requires  the  total  traction  to  vanish  there,  i.e. 


nt=  0  ,at  x^  =  f{x^,x^,t)  .  (11) 

The  above  condition  also  applies  to  a  continuous  material.  Reportedly  (see  Johnson  and 
Jackson  1987)  when  this  boundary  condition  is  applied  to  problems  concerning  granular 
media,  it  tends  to  introduce  a  singular  behavior  of  the  traction  in  the  neighborhood  of 
the  interface  that  is  physically  unrealistic,  as  the  stress  tensor  vanishes  rapidly  through  a 
thickness  that  is  typically  less  than  the  thickness  of  one  particle  diameter.  However,  obtaining 
“depth-averaged”  solutions  is  intended  ,and  this  difficulty  should  not  be  a  problem,  except 
at  the  flow  boundaries. 
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Furthermore,  it  is  important  to  discuss  the  meaning  of  a  “free  surface”  in  a  physical 
sense.  As  discussed,  it  is  generally  accepted  that  an  avalanche  may  he  considered  to  be  a 
flow  consisting  of  three  distinct  regimes.  The  dense  core  of  the  flow  is  normally  obscured 
by  a  “dust  cloud,”  which  is  a  layer  of  highly  suspended  particulate.  From  observation,  this 
layer  is  dominated  by  the  interstitial  fluid  behavior  (Lang  et  al.  1989).  Intuitively,  if  one 
defines  the  “free  surface”  as  the  interface  of  the  core  of  the  flow  with  the  turbulent  dust 
cloud,  it  would  be  more  appropriate  to  consider  the  drag  at  this  rough  interface  rat!ier  than 
a  vanishing  value  of  traction.  However,  it  can  be  demonstrated  (McClung  and  Schaerer  1983) 
that  inclusion  of  such  a  term  has  a  negligible  effect  on  the  overall  resistance  to  motion  of  the 
mass. 

With  these  considerations  in  mind,  the  stress-free  condition  is  assumed,  and  the  resulting 
boundary  conditions  at  the  free  surface  are 


written  in  the  and  x^  directions  respectively. 

At  the  base  of  the  flow,  the  surface  may  be  described  by  a  function  b  =  b{x^  ,  x^),  such 
that  the  function  b  may  be  arbitrarily  definea  mr  a  given  avalanche  bed. 

The  tangency  of  flow  at  the  base  is  expressed  by  the  kinematic  condition 

un  =  0  ,at  x^  =  6(x^,x^)  .  (13) 

Expanding  the  above  in  three  spatial  dimensions  yields 

db  f  dC  \  \  db  1  0 

^v<i>  +  (^1  +  ^/j  ^^^<2>  -  v<3>J  =  0  x^  =  6(x\x^)  .  (14) 

Applying  the  assumption  that  the  basal  boundary  layer  is  of  a  thickness  proportional 
to  only  a  few  particle  diameters,  this  layer  is  neglected.  This  cissumption  greatly  simplifies 
the  mathematical  formulation,  as  in  this  layer  collisional  and  frictional  stresses  must  be 
considered.  An  appropriate  boundary  condition  has  been  derived  by  Hui  et  al.  (1984)  for 
such  a  case,  based  on  kinetic  theory.  However,  some  knowledge  of  the  specularity  of  particle- 
wall  collisions  is  required. 

By  considering  only  cohesionless,  dry  granulate,  the  basal  shear  traction  may  be  expressed 
as 


nt  —  n{nt  •  n)  = -nt- n- — j-tan6  — D„lul  v  ,at  x^  =  6(x^,x^)  . 

I  U.  \ 


where  the  first  term  on  the  right  hand  side  is  the  Coulomb  sliding  law,  and  tan  6  is  the  bed 
friction  angle  whose  value(s)  depend  on  the  nature  of  the  particles,  degree  of  packing,  and 
the  contact  surface.  The  remaining  term  on  the  right-hand  side  is  the  boundary  drag  term. 
It  should  be  noted  that  the  above  is  compatible  wuth  the  kinematic  condition  at  the  base  (eq 
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13).  When  the  basal  surface  agrees  with  the  coordinate  generating  surface,  n  =  and 
the  kinematic  condition  at  the  base  eq  13  reduces  to 

w<3>  =  0  .  (16) 

Equation  15  at  the  base  reduces  to 


*<3  1> 


ftan^ 

iT7r*<3'>+ 


D,\ 


(17) 


^<3  2>  =  -V<2> 


tan  6 


*<3  3>  +  I  V 


In  view  of  the  experimental  work  (Lang  et  al.  1989),  which  encompassed  only  the  case 
of  the  basal  surface  replicating  the  coordinate  surface,  the  following  discussion  will  focus  on 
this  limiting  case.  Also,  the  scaling  arguments  will  demonstrate  the  inconsequence  of  weak 
curvature.  For  strong  curvature  in  the  lateral  direction,  the  theory  may  be  easily  reduced  to 
channel  flow  equations. 


Limiting  cases  of  the  governing  equations 

The  governing  equations  that  were  presented  are  now  nondimensioneilized  and  averaged 
over  the  depth  and  width,  in  order  to  provide  a  set  of  tractable  equations.  The  nondi- 
mensionalization  scheme  allows  for  any  number  of  characteristic  length  scales,  that  differs 
from  more  conventional  approaches.  Becker  (1976)  presents  a  thorough  treatment  of  this 
method,  commonly  termed  a  “distorted  model.”  However  he  prefers  the  more  distinctive 
nomenclature  of  “Pythagorean  method,”  as  the  Pythagoreans  were  motivated  by  the  idea  of 
geometrical  proportion  being  the  basis  of  all  things.  By  applying  this  method,  inspectional 
analysis  indicates  limiting  forms  of  the  equations  that  may  be  more  easily  treated.  The  ap¬ 
plication  of  this  method  has  been  accepted  as  the  appropriate  procedure  when  attempting  to 
model  large-scale,  geophysical  phenomena.  By  adopting  the  concept  of  a  “distorted  model,” 
consideration  of  more  than  one  characteristic  length  scale  is  edlowable.  Such  an  approach  is 
commonly  used  and  applied  to  large-scale,  free-surface  flow  phenomena  such  as  atmospheric 
and  oceanic  motions.  A  specific  example  is  the  well-known  “shallow  water  model.”  Scaling 
provides  the  means  to  neglect  small  variations  in  topography  that  do  not  significantly  affect 
the  lateral  or  longitudinal  spread  of  the  flowing  material.  Curvature  in  the  lateral  direction 
will  be  neglected,  as  discussed  in  a  later  section.  For  strong  curvature  in  the  lateral  direction, 
the  equations  are  easily  reduced  to  a  channel  flow  model. 

Also,  assumptions  are  made  to  simplify  the  balance  equations  and  boundary  conditions. 
Many  of  these  result  directly  from  the  scaling  arguments,  and  others  are  based  on  experience 
with  the  physical  problem  at  hand  The  nondimensionalized  equations  are  then  averaged 
over  the  depth  and  width  of  the  flow.  This  is  a  classical  approximation  method,  commonly 
applied  in  evaluating  boundary  layer  flows.  The  result  is  a  set  of  equations  that  represent  a 
balance  between  the  average  inertial  and  resistant  forces  for  each  position  in  the  longitudinal 
direction.  Ordering  of  arguments  allows  further  simplification.  The  need  to  consider  the 
detailed  character  of  the  flow  is  thus  avoided,  and  sufficient  information  for  an  engineering 
application  is  obtained. 
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Nondimensionalization  scheme 

At  this  point,  a  scheme  for  nondimensionalization  of  the  governing  equations  is  presented. 
This  step  of  modeling  geophysical  phenomena  is  always  to  some  extent  a  “leap  of  faith.” 
That  is  to  say,  choosing  an  appropriate  model  depends  on  experience,  intuition  and  explicit 
dynamic  reasoning.  The  essential  dynamic  character  of  the  geophysical  system  must  be 
retained  but  the  validity  of  the  model  is  examined  after  the  fact. 

Although  the  depth,  length  and  width  of  an  avalanche  varies  spatially  and  temporally,  it 
is  assumed  that  characteristic  values  O  may  be  sensibly  chosen  such  that 


x'  =  Cx\  (18) 

where  the  x'  are  the  nondimensional  lengths.  The  concept  of  an  avalanche  as  a  gravity 
current  suggests  the  following  representation. 


C'  r  iii/2 

^<t>  =  ^  J  =  V’w<i>  .  (19) 

The  above  may  be  interpreted  with  respect  to  the  fluid  dynamic  concept  of  “head”  or  the 
free  fall  of  a  mass  point  in  a  gravity  fleld  through  a  given  distance,  although  the  length  scale 
may  differ.  This  leads  to  an  expression  for  time  as 


t  =  Ti  = 


(20) 


As  the  Proude  number  is  the  governing  parameter  for  the  problem  under  consideration,  the 
stresses  are  scaled  according  to  a  hydrostatic  relationship, 


t<t  i>  —  PoQ^  COS^o^<»»>  (21) 

t<tj>  =  poQC,  sin^o?<»j>>  (22) 

where  the  null  subscript  denotes  a  characteristic  value  for  the  quantity  under  consideration, 
while 


sin  Co  =  cos  tan  5  .  (23) 

The  expression  for  sin  Co  is  a  logical  choice,  as  both  the  driving  force  and  the  resistance  are 
represented  in  eq  23,  which  both  influence  the  order  of  magnitude  of  the  resultant  shear 
stresses.  The  angle  C  may  also  be  characterized  as 


O  - 

(24) 

Also,  the  following  definitions  are  applied; 

£* 

(25) 

C} 

^  =  ^  ■ 
n 

(26) 

In  order  to  evaluate  the  basal  boundary  conditions  the  drag  coefficient  is  defined  in 
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terms  of  the  Voellmy  coeflScient  w  ,  which  is  traditionally  applied  to  flowing  snow,  as  follows; 

=  —  .  (27) 

G7 

w  has  the  units  of  m/s^  such  that  it  is  a  natural  choice  to  scale  it  with  the  gravitational 
constant  such  that  w  is  related  to  the  dimensionless  form  w  by 


G7  =  Qw 


(28) 


Depth  and  width  averaging 

When  a  functional  form  of  solution  to  a  set  of  governing  equations  cannot  be  found  (i.e. 
a  solution  that  satisfles  the  equations  identically  at  each  point  in  space  and  time  and  that 
approaches  the  appropriate  boundary  values)  it  may  be  advantageous  to  try  and  satisfy  the 
equations  on  the  average,  if  the  geometry  of  the  problem  allows  for  this  approach  to  result  in 
a  meaningful  solution.  If  the  equations  may  be  integrated  with  respect  to  across  the  depth 
of  flow  and  with  respect  to  across  the  width  of  flow,  the  resulting  equations  will  represent 
a  balance  between  the  average  inertial  forces  and  resistant  forces  for  each  location.  In 
this  manner,  a  velocity  distribution  may  be  obtained  that  satisfies  this  averaged  balance  of 
forces  but  that  does  not  satisfy  the  balance  at  each  point  across  the  height  and  width  of  the 
flowing  mass. 

This  method  is  a  classical  approximation  method  that  is  widely  used  in  examining  bound¬ 
ary  layer  behavior,  and  was  first  introduced  by  von  Karmin  (1921)  and  was  further  improved 
by  Pohlhausen  (1921).  The  results  that  are  obtained  from  such  an  approximation  are  found 
to  be  reasonably  accurate  in  many  instances.  By  acknowledging  thatv<i>  and  v<2>  are  con¬ 
tinuous  on  the  interval  [  0,  /],  the  first  mean  value  theorem  of  calculus  ascertains  the  existence 
of  a  point  p  contained  in  the  open  interval  (0,  / )  such  that 


Jo 


®<i>  dx 


=3 


=  t'<.>(p)  /  dx^ 

JQ 


=  hv 


<»>  > 


(29) 


and  furthermore,  it  may  be  assumed  that 


(30) 


where  the  notation  »<<>  implies  the  “depth-averaged”  value  of  «<,•>  and  f  —  0  =  h,  the 
nondimensional  flow  height.  The  latter  equation  is  exact  if  the  velocity  profile  over  the  depth 
of  the  flow  is  uniform. 
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Let 


€3  i€3  2  -^2  3  =  0(^3  l)  (31) 

and  retain  only  terms  of  0(€3  1)  .  The  depth-averaged  velocity  squared  is 

^  y  (32) 

where  the  Boussinesq  momentum  coefficient  will  be  taken  to  be  near  unity  such  that 


v^<i>  ^  •  (33) 

To  be  realistic,  the  Boussinesq  coefficient  is  greater  than  one,  but  this  approximation  is 
sufficient  when  the  velocity  profile  is  blunt.  For  a  power  law  velocity  distribution  (see, 
for  example.  Burgers  1924)  that  typifies  a  pure  Newtonian  turbulent  fiuid,  the  Boussinesq 
momentum  coefficient  is  1.05  (Savage  1989b).  The  normal  stress  at  the  bed  is  then 


i<3  3> 


cos  Co 


dC 


(34) 


fS=0 

The  use  of  an  “earth  pressure  coefficient”  is  employed,  such  that  the  normal  stresses 
parallel  and  perpendicular  to  the  path  may  be  equated  by  employing  an  euth  pressure 
coefificient  koctpajj 


f<l  1>  —  f^actpaa$i<3  3>  •  (35) 

kactpaaa  has  a  value  of  kact  or  kpassi  determined  by  the  spatial  evolution  of  the  pile  as  either 
elongation  (active  case)  or  contraction  (passive  case)  with  respect  to  the  flow  direction  x^. 
From  a  standard  Mohr’s  circle  interpretation,  it  may  be  shown  that 

1  T  ^1  -  (1  tan^  6)  cos^  (^1  /cos^  ^  -  1 ,  (36) 


It  is  also  assumed  that 


and 


for  dv<i>/da}  = 


_  j  fh 

V<l>V<2>  =  jj  V<l>V<2>dx^ 


V<l>V<2>  -  5<i>W<2>- 


Also 


(37) 


(38) 


23  2 


(39) 
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The  integrand  f<2  i>  is  neglected  by  ordering  A2  i  €32  «  0(€i  2),  as  this  shear  stress 
would  only  be  locally  effective  in  the  case  of  confined  narrow  channel  flows.  It  has  been  shown 
(Savage  1989  and  Lang  1988)  that  sidewall  shearing  may  be  incorporated,  if  necessary,  as  an 
effective  friction. 

Also 

_  1  _ 

v<i>v<2>  =  ^  «<i>v<2>dx^  (40) 

and  assume  that 


v<i>w<2>  = 


By  neglecting  the  2>  integrand,  the  lateral  pressure  gradient  may  be  integrated  as 


d  d  ' 


^actpaas^^  /  dC 

- - ;r-  1  cosC  +  IC-^v 

2  cos  Co  \  dx^ 


]] 


based  on  the  value  of  the  t<3  3>  stress  component.  An  approximation  for  the  t<3  2>  stress 
may  be  made  from  eq  42  from  relating  the  shear  to  normal  stress  such  that 


L 


v<2>h'^ 

I  A3  2  cos  Co  tan  6 


The  nondimensionalized  depth-averaged  equations  now  may  be  further  simplified  by 
width  averaging.  By  defining  the  (nondimensional)  half-width  of  the  pile  as 

bm  =  bm{x^,i)  (44) 

and  assuming  that  the  pile  of  granulate  retains  a  parabolic  shape  in  the  x^  plane,  the  height 
may  be  expressed  as 


h{x^,  x^,  t)  =  hoix^,  t) 


(45) 


where  /io(  t )  is  the  centerline  depth  (refer  to  Figure  1).  This  assumption  is  fully  supported 
by  the  experimental  results  discussed  in  Lang  et  al.,  1989.  Although  natural  avalanches  do 
not  retain  such  perfect  symmetry  due  to  uneven  topography  and  particle  size  variations, 
some  degree  of  symmetry  is  usually  observed  in  natural  avalanche  deposition  and  a  parabolic 
distribution  is  reasonably  replicated. 


1  I  r 


-bm(x^,t)  0  6m(x^^) 

Figure  1:  Approximation  of  pile  shape  in  the  x^  plane 
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Figure  2:  Lateral  pile  velocity 


Assuming  parabolic  symmetry  of  the  piie  implies  that  v<2>  ,  the  nondimensional  depth 
averaged  lateral  pile  velocity,  varies  linearly  with  respect  to  as  in  Figure  2  such  that 


V<2>  = 


dbm 
dt  brn 


(46) 


Integrating  the  continuity  eq  3  over  the  depth  and  half-width,  applying  Leibniz’  rule,  and 
invoking  the  symmetry  conditions  in  (45  and  46)  gives 


9  T=  .2  ^  /'*’”*!  j-2  d(^  d  /‘*™  k‘‘  ,_2  - 


d 

dx^  Jo 

Define  a  “width” -averaged  depth  h 


dC  d  .r 


1  /'*"•  - 

h  =  ~  h  dx^ 
Om  Jo 


(47) 


(48) 


By  assuming  that  the  depth  has  a  parabolic  distribution  in  the  =  constant  planes 

h  =  .  (49) 

Furthermore,  a  depth-  and  width-averaged  velocity  may  be  defined  again  by  use  of  the  mean 
value  theorem  of  calculus  such  that 


1  /■*  -  2 
v<i>  =  TT  /  v<i>hdx^  . 

Omft  Jo 


With  these,  continuity  becomes 
d 


9  {hobjn  v<i>)  ,  9(ho6m),2  dQ  d  ( h*  l>m )  _  « 


(50) 


(51) 


Width-averaging  of  the  momentum  balance  equations  requires  an  approximation  for  the 
norm  of  the  velocity  vector.  By  acknowledging  that  sidewise  spreading  contributes  in  a 
neglible  and  fairly  constant  maimer  to  the  value  of  the  depth-averaged  velocity,  one  approx¬ 
imation  to  oe  considered  is 


Furthermore,  let 


Hi-  ®<i> 


(52) 
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A  = 
B  = 

C  = 


€32  kactpass  I 


fC 


dx^ 


K,  €3  1 


dx^ 


Fi  = 
Fi  = 


€3 1  ro 
1 

€32  ^ 


/  hv<i>v<2>  dx^  —  v<i>  f  hv<^2>  dx^ 

Jo  Jo 


(53) 


where  the  last  approximation  is  physically  realistic,  as  the  longitudinal  velocity  has  been 
assumed  to  have  a  blunt  profile  over  depth  and  does  not  appear  to  vary  over  the  width  at 
any  one  transverse  section.  Applying  the  aforementioned  definitions  and  approximations  to 
the  momentum  equation  in  the  e<2>  direction  eq  6,  and  by  application  of  Leibniz’  rule, 
termwise  integration  from  x^  =  0  to  =  h  and  =  0  to  =  bm  gives 


+  (cosC  +  Bfi<i>) -^^(«<2>  I  )  ^l-^)tan^(cosC  +  5fi|i>) 

6m 

-  -^^v<i>(v<2>  j  )  ^1 - (54) 


€2  1 


brnho 


X  ^/  = 

)  +  ^(^<2> 


,  .  Cho  d 


bm  bm  bm 

The  above  equation  is  easily  obtained  by  simple  integration  by  noting  that 


15  at  ^  ■ 


h 


=  0 


V<2>  j 

o 


0  . 


Similarly,  let 


G  — €3  1  kactpaas 


(55) 


By  the  same  method  the  width-averaged  momentum  balance  in  the  e<i>  direction  is 
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3G  d  26mftoCOsC  8bmhlBv\y^  / 

Abmhodx^  L  ^  15  J  r 


AhoC 


^  j^cos Ctan(5  +  5i;^j^  tan 5 j 


Fiv^ 


i> 


1  ^brnho 


2hoC\ 

5  ) 


sinC 


_  dv<i>  .  dv^i>  2Cho  dv<i> 

dt  Q^i  5  Of 

^  dv<i> 
di 


(56) 


Equations  51,  56,  and  54  comprise  a  complete  boundary  value  problem  for  the  free  surface 
flow  of  a  cohesionless  granulate  down  an  open  terrain  of  variable  curvature  in  the  longitudinal 
direction.  The  requirement  of  pile  symmetry  has  been  imposed,  which  is  satisfled  by  the 
experiments  in  Lang  et  al.  (1989)  that  meet  the  criteria  of  being  friction-governed  flows. 


NUMERICAL  SOLUTIONS  AND  COMPARISON  TO  EXPERIMENTS 

Because  a  closed  form  solution  is  not  available  for  the  system  of  equations  which  was 
outlined,  a  numerical  solution  must  be  sought.  Savage  and  Hutter  (1989)  show  that  a  La- 
grangian  formulation  is  the  natural  choice  for  this  system  of  equations  in  order  to  determine 
the  position  of  the  granulate  surface  through  time.  With  this  method,  the  boundaries  of 
the  finite  difference  grid  are  convected  with  the  depth-  and  width-averaged  velocities  of  the 
granular  material. 

Two  different  models  are  developed.  First,  the  boundary  drag  term  is  neglected  and 
predictions  are  based  solely  on  modelling  the  material  with  the  Coulomb  friction  law.  Then, 
it  is  demonstrated  that  a  flow  transition  regime  exists  where  the  boundary  drag  term  must 
be  included  in  the  model  in  order  to  predict  the  position  of  the  pile  through  time,  and  more 
specifically,  the  final  runout  zone  position.  An  artificial  viscosity  term  fid^u/dx^  is  added  to 
the  equation  of  motion  in  order  to  smooth  the  solution. 


Unconstrained  Coulomb  flow  model  with  constant  bed  friction 

Now  the  depth-  and  width-averaged  balance  equations  are  solved.  In  addition  to  pre¬ 
dicting  the  length,  height  and  longitudinal  velocity  through  time,  the  sidewise  spreading  and 
rate  of  sidewise  spreading  may  be  determined.  The  results  from  channel  flow  models  (Lang 
1988)  are  used  as  a  guideline.  A  first  approximation  is  applied,  where  it  is  assumed  that 
the  boundary  drag  term  may  be  neglected,  and  that  d  remains  constant  through  the  length 
of  the  pile.  This  provides  the  simplest  set  of  width-averaged  equations.  For  this  system  of 
equations,  the  condition  is  applied  that  C  <  <  1.  The  depth-  and  width-averaged  equations 
as  presented  are  then  simplified.  For  comparison  with  simple  laboratory  findings  a  constant 
bed  friction  angle  provides  good  results.  This  model  was  tested  against  the  experimental 
results  presented  in  Lang  et  al  1989.  The  numerical  model  and  laboratory  results  show  rea¬ 
sonable  agreement.  For  example,  experiment  DA029  was  released  at  an  initiation  zone  angle 
of  35° .  DA029  consisted  of  20.0  kg  of  2-  to  3-  mm  quartz  particles,  reaching  a  maximum 
velocity  of  3.6  m/s  at  0.9  s.  Figure  3  shows  the  experimental  values  of  velocity  and  position 
for  the  leading  and  trailing  edges  through  time  as  compared  with  the  theoretical  predictions. 
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Figure  3:  Theoretical  predictions  and  experimental  values  of  nondimensional  position  and 
velocity  of  the  leading  and  trailing  edge  as  a  function  of  nondimensional  time  for  experiment 
DA  029 

The  leading  edge  is  the  point  on  the  =  0  axis  that  travels  the  farthest,  and  the  trailing 
edge  is  the  point  on  the  =  0  axis  that  travels  the  least  amount  of  distance  through  time. 
Due  to  the  method  of  photography  (see  Lang  et  al  1989),  a  nondimensional  time  error  of 
0.52  may  be  incurred  in  determining  the  initial  frame,  which  may  result  in  some  of  the  offset. 
Nevertheless,  the  theory  predicts  the  general  characteristics  of  the  longitudinal  progress  of 
the  granular  mass,  and  the  final  position  of  the  forward  edge  in  the  runout  zone  is  accurately 
predicted.  The  prediction  of  the  position  and  velocity  of  the  leading  edge  of  an  avalanche 
is  the  critical  information  required  to  calculate  the  impact  pressure  on  a  structure  and  to 
motivate  zoning  restrictions. 

In  Figure  4  the  progression  of  the  width  and  height  at  selected  time  increments  is  pre¬ 
sented.  The  time  increments  were  chosen  to  correspond  with  frame  times.  The  left  hand  side 
of  the  figure  represents  the  progression  of  the  half  width  of  the  model  and  experiment  (see 
Lang  et  al  1989)  at  selected  time  increments,  and  the  right  hand  side  is  the  predicted  depth 
profile  through  time.  Due  to  limited  space  in  the  laboratory,  it  was  not  possible  to  record 
the  height  of  the  pile  through  time,  but  the  general  features  are  reasonably  reproduced  by 
the  theory;  that  is,  the  pile  collapses  rapidly,  elongates  and  flattens  to  a  maximum  prior  to 
encountering  the  runout  zone.  Then  as  the  forward  portion  of  the  pile  stops,  the  trailing 
edge  retains  momentum  and  piles  up  over  the  static  portion.  When  the  angle  of  repose  is 
exceeded,  the  pile  relaxes  back  to  the  internal  angle  of  friction,  and  the  tail  comes  to  rest 
at  some  angle  less  than  the  angle  of  repose  (Fig.  4).  The  maximum  final  (nondimensional) 
pile  height  was  measured  to  be  0.21,  and  the  model  predicts  0.23.  The  model  overestimates 
the  (dimensional)  lateral  spread  of  the  pile,  by  approximately  9%,  and  the  overestimation  of 
the  lateral  spread  causes  the  underestimation  of  the  pile  length  through  time,  particularly  in 
the  runout  zone.  The  pile  is  predicted  to  spread  more  rapidly  in  the  lateral  direction  in  the 
initiation  zone. 

For  experiments  that  were  conducted  at  an  initial  bed  angle  (or  starting  zone  angle)  of 
35°,  the  constant  bed  friction  angle  model  is  a  good  predictor.  For  laboratory  simulations 
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Figure  4:  Sequence  of  nondimenstonal  computed  and  experimental  pile  width  and  computed 
pile  height  at  frame  time  increments  for  DA029 
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conducted  with  greater  starting  zone  angles,  the  model  greatly  overpredicts  the  final  position 
of  the  granulate  mass  in  the  runout  zone.  For  example,  experiments  DA025  and  DA033  were 
also  modeled  with  the  constant  bed  friction  model.  Both  of  these  experiments  were  released 
at  an  initiation  angle  of  45°.  For  DA025,  20.0  kg  of  2-  to  3-  mm  quartz  particles  were  used. 
The  resulting  predictions  of  position  and  velocity  through  time  show  that  the  final  positions 
of  the  leading  edge  and  the  tail  are  greatly  over  predicted.  Position  and  velocity  profiles  for 
DA033  18.4  kg  of  2-3  mm  marmor  show  also  that  the  toe  and  tail  positions  are  overpredicted 
with  this  model  (for  details,  see  Lang  1992). 

This  result  agrees  with  problems  encountered  when  attempting  to  model  natural  snow 
avalanches.  Gubler  (1987)  discusses  the  problem  of  modeling  higher  speed  avalanches  with 
traditional  avalanche  models.  The  friction  coefficients  must  be  increased  by  a  factor  of  2,  in 
general,  in  order  to  force  the  avalanche  to  come  to  rest  in  the  appropriate  position  in  the 
runout  zone. 

Because  this  numerical  model  accurately  predicts  the  experiments  conducted  at  35°  which 
results  in  lower  velocities  and  accelerations  than  for  larger  angles,  it  is  possible  that  a  tran¬ 
sition  in  flow  regime  occurs  in  the  boundary  layer  where  the  boundary  drag  force  becomes 
non-negligible  at  higher  velocities.  Its  inclusion  in  the  numerical  model  is  presented  in  the 
following  section.  Gubler  (1987)  discusses  a  similar  flow  regime  transition  for  natural  snow 
avalanches  from  flexible  sliding  body  to  partial  fluidization.  Furthermore,  Buser  and  Frutiger 
(1980)  found  that  the  inclusion  of  both  sliding  bed  friction  and  boundary  drag  underpredict 
velocities  for  very  large  natural  snow  avalanches.  This  result  may  also  be  an  indication  of 
flow  transition  behavior,  in  the  sense  that  an  avalanche  may  initiate  in  the  flexible  sliding 
body  type  flow  and  reach  a  transitive  velocity  where  partial  fluidization  occurs.  If  the  ad¬ 
ditional  resistance  due  to  boundary  drag  is  applied  throughout  the  flow  time,  the  theory 
would  underpredict  velocities.  However,  if  a  critical  velocity  could  be  determined  that  would 
predict  the  transition,  and  hence  the  appropriate  application  of  an  additional  boundary  drag 
term,  perhaps  this  problem  could  be  circumvented. 


Unconstrained  Coulomb  flow  with  boundary  drag  term 


It  was  mentioned  in  the  previous  discussion  that  the  constant  bed  friction  angle  model 
becomes  a  poor  predictor  for  initiation  zone  angles  greater  that  35° .  At  this  point  it  seems 
sensible  to  attempt  to  model  the  additional  resistance  to  flow  by  inclusion  of  the  boundary 
drag  term.  Again,  average  values  are  used  to  avoid  the  possibility  of  the  pile  height  approach¬ 
ing  zero.  Also,  the  depth-  and  width-averaged  drag  term  includes  the  width  of  the  pile  in  the 
denominator,  that  is  zero  at  the  midpile.  Therefore,  the  following  average  of  the  drag  term 
is  included: 


2  €3 1  w  hav  bav 


(57) 


As  an  example,  experiment  DA025  was  modeled  with  the  inclusion  of  the  drag  term,  and 
results  are  presented  in  Figures  5  and  6.  DA025  was  conducted  at  an  initiation  angle  of  45° 
and  20.0  kg  of  2-3  mm  quartz  particles  were  used.  Figure  5  depicts  the  resulting  predictions 
of  position  and  velocity  through  time.  The  final  position  of  the  leading  edge  is  somewhat 
underpredicted,  but  the  general  motion  of  the  pile  through  time  is  well  mirrored.  The  model 
does  not  allow  the  maximum  leading  edge  velocity  to  be  attained.  It  may  be  restated  that 
ideally,  the  boundary  drag  term  should  perhaps  be  applied  only  when  a  critical  velocity  (or 
critical  Froude  number)  is  attained. 
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Figure  5;  Theoretical  predictions  and  experimental  values  of  nondimensional  position  and 
velocity  of  the  leading  and  trailing  edge  as  a  function  of  nondimensional  time  for  experiment 
DA  025 

As  shown  in  Figure  5,  the  initial  spread  of  the  pile  is  overpredicted,  but  in  the  runout  zone, 
the  shape  is  quite  accurately  replicated.  The  pile  is  predicted  to  be  shorter  in  the  longitudinal 
direction  and  slightly  wider  than  the  actual  pile  shape.  The  final  maximum  (nondimensional) 
pile  height  is  extremely  accurate:  the  theory  predicted  0.205  and  the  experimental  result  was 
0.206. 

CONCLUSIONS  AND  REMARKS 

Governing  equations  and  closure  relations  were  proposed  to  represent  the  flow  of  a  co¬ 
hesionless  granular  medium  on  any  given  topography.  The  goal  was  to  develop  an  easily 
applied  model  to  predict  the  bulk  flow  characteristics  of  the  mass  of  granulate.  Indeed,  the 
anticipated  result  of  every  avalanche  model  is  to  predict  the  final  runout  zone  position  of 
the  mass  and  the  possible  impact  pressures  on  structures  in  the  path  of  an  avalanche.  A 
purely  mechanical  theory  was  proposed  and  the  granular  mass  wais  treated  as  a  continuum. 
The  balance  equations  in  terms  of  physical  components  in  an  adaptable  and  appropriate 
coordinate  system  were  derived.  A  constitutive  assumption  was  justified  by  considering  the 
different  flow  regimes  in  question,  and  the  experimental  results  from  Lang  et  al.  (1989),  and 
other  investigations. 

The  governing  equations  were  nondimensionalized  and  averaged  over  the  depth  and  width, 
in  order  to  provide  a  set  of  tractable  equations.  The  nondimensionalization  scheme  allows  for 
any  number  of  characteristic  length  scales,  which  differs  from  more  conventional  approaches. 
By  applying  this  method,  inspectional  analysis  indicates  limiting  forms  of  the  equations 
that  may  be  more  easily  treated.  The  application  of  this  method  has  been  accepted  as 
the  appropriate  procedure  when  attempting  to  model  large  scale,  geophysical  phenomena. 
Assumptions  were  made  to  simplify  the  balance  equations  and  boundary  conditions,  which 
resulted  directly  from  the  scaling  arguments  or  were  based  on  experience  with  the  physical 
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Figure  6;  Sequence  of  nondimensional  computed  and  experimental  pile  width  and  computed  pile  height  at  frame 
time  increments  for  DA025 
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problem  at  hand.  Changes  in  the  volume  fraction  of  the  granulate  and  variations  in  topog¬ 
raphy  in  the  lateral  direction  wer»;  neglected.  The  nondimensionalized  equations  were  then 
averaged  over  the  depth  and  width  of  the  flow,  a  classical  approximation  method,  commonly 
applied  in  evaluating  boundary  layer  flows  resulting  in  a  set  of  equations  which  represent  a 
balance  between  the  average  inertial  and  resistant  forces  for  each  position  in  the  longitudi¬ 
nal  direction.  With  this  method,  no  information  can  be  provided  to  determine  the  detailed 
character  of  the  flow,  but  sufficient  information  for  an  engineering  application  is  obtained. 

A  Lagrangian  discretization  scheme  was  applied  in  order  to  numerically  evaluate  the 
magnitude  of  velocity,  and  predict  the  spatial  evolution  of  the  granular  pile  through  time. 
This  type  of  analysis  allows  for  the  particulate  to  be  advected  with  the  finite  difference  grid, 
and  was  shown  by  Savage  and  Hutter  (1989)  to  be  the  appropriate  choice  for  such  a  problem. 

Two  models  were  developed  for  the  free  surface  flows.  The  first  assumed  a  constant  bed 
friction  angle,  and  neglected  any  boundary  drag  effects.  This  model  proved  to  be  adequate 
in  predicting  the  experiments  discussed  in  Lang  et  al.  (1989),  which  were  conducted  at  bed 
surface  angles  of  35",  but  failed  for  larger  values  of  initiation  zone  angles.  The  method  of 
photography  incurred  a  nondimensional  time  error  of  0.52  in  determining  the  initial  frame, 
which  may  have  resulted  in  some  of  the  inaccuracy.  Nevertheless,  the  theory  accurately 
predicts  the  general  characteristics  of  the  longitudinal  progress  of  the  granular  mass,  the 
velocity  with  time,  and  the  final  position  of  the  forward  edge  in  the  runout  zone  for  the 
35"  experiments  (e.g.,  DA029).  The  progression  of  the  height  and  width  at  selected  time 
increments,  corresponding  with  frame  times,  was  also  well  reproduced  by  the  theory.  The 
pile  collapses  rapidly,  elongates  and  flattens  to  a  maximum  prior  to  encountering  the  runout 
zone.  Then  as  the  forward  portion  of  the  pile  steps,  the  trailing  edge  retains  momentum 
and  piles  up  and  over  the  static  portion.  When  the  angle  of  repose  is  exceeded  (Fig.  4) 
the  tail  reverts  back  to  an  angle  equivalent  to  or  less  than  the  internal  angle  of  friction. 
The  final  angle  (actual  or  predicted)  depends  on  the  amount  of  backwashing  that  occurs. 
The  maucimum  final  (nondimensional)  pile  height  was  predicted  to  within  8%.  The  model 
overestimates  the  lateral  spread  of  the  pile  and  this  causes  the  underestimation  of  the  pile 
length  through  time,  particularly  in  the  runout  zone.  The  pile  is  predicted  to  spread  more 
rapidly  in  the  initiation  zone.  But  the  average  flow  procession  of  the  pile  is  well  predicted. 

For  laboratory  simulations  conducted  with  greater  starting  zone  angles,  the  model  that 
neglects  any  flow  resistance  due  to  drag  overpredicts  the  final  position  of  the  granulate  mass 
in  the  runout  zone.  This  result  agrees  with  problems  encountered  when  attempting  to  model 
natural  snow  avalanches.  Attempting  to  model  higher  speed  avalanches  with  a  basic,  or  all- 
encompassing,  avalanche  model,  and  producing  erroneous  results  is  not  a  new  or  surprising 
problem.  Gubler  (1987)  discusses  the  problem  of  incre:<  >ing  friction  coefficients  by  a  factor  of 
2,  in  general,  in  order  to  force  avalanche  models  to  come  to  rest  in  the  appropriate  position 
in  the  run-out  zone,  when  attempting  to  model  snow  avalanches  that  attain  higher  velocities. 

Since  the  first  model  breaks  down  for  initiation  zone  angles  greater  that  35"  ,  an  attempt 
was  made  to  model  the  additional  resistance  to  flow  by  inclusion  of  the  boundary  drag  term. 
Experiment  DA025  was  modeled  with  the  inclusion  of  the  boundary  drag  term,  and  results 
were  presented  in  Figures  5  and  6.  This  experiment  was  released  at  an  initiation  angle  of 
45". 

For  DA025,  the  final  position  of  the  leading  edge  was  underpredicted,  but  the  general 
motion  of  the  pile  through  time  is  well  mirrored.  The  model  does  not  allow  the  maximum 
leading  edge  velocity  to  be  attained.  It  is  not  an  unreasonable  idea  that  the  boundary  drag 
term  should  perhaps  be  applied  only  when  a  critical  velocity  (or  critical  Froude  number)  is 
attained.  Also,  the  initial  spread  of  the  pile  was  overpredicted,  but  in  the  runout  zone,  the 
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shape  is  quite  accurately  replicated.  The  pile  is  predicted  to  be  shorter  in  the  longitudinal 
direction  and  slightly  wider  than  the  actual  pile  shape.  The  final  maximum  (nondimen- 
sional)  pile  height  was  accurately  predicted,  with  a  theoretical  value  of  0.205  compared  to 
the  experimental  result  of  0.206. 

In  general,  the  numerical  models  and  experiments  show  good  agreement  and  provide 
valuable  information.  These  results  are  justification  for  the  characterization  of  the  flow,  and 
persuance  of  averaged  quantities.  The  resulting  predictions  of  geometry  and  velocity  through 
time  are  encouraging.  However,  some  subjectivity  is  still  required  to  determine  appropriate 
constitutive  parameters  and  it  is  unknown  if  the  models  can  represent  naturally  occurring 
events. 

Future  work  should  include  testing  these  models  against  naturally  occurring  events,  in 
particular,  snow  avalanches,  since  experimental  data  are  available  (Gubler,  1987).  This  would 
be  no  minor  task,  considering  one  would  have  to  attempt  to  determine  in-situ  properties  of 
the  snowpack.  Also,  the  model  presented  was  developed  to  include  curvature  variations  of  the 
slide  path  in  the  lateral  direction.  If  one  were  to  completely  test  the  model  against  a  naturally 
occurring  event,  it  would  have  to  be  determined  how  important  the  topographical  variations 
in  the  lateral  direction  are  in  afiFecting  the  final  runout  position.  Certainly,  constraining  the 
flow  in  a  channel  prohibits  dilatation  to  some  degree.  This  would  entail  more  experimental 
work  on  simple  materials.  It  may  be  sufficient  to  characterize  a  given  natural  topography 
by  simply  combining  the  channel  flow  and  unconstrained  flow  models,  and  specifying  the 
channel  width  as  a  function  of  longitudinal  direction,  and  the  position(s)  where  channel  flow 
is  no  longer  appropriate.  Asymmetry  in  the  spreading  direction  may  or  may  not  be  a  decisive 
factor  in  determining  the  final  runout  position  of  the  flow. 

Also,  in  order  to  compare  the  model  to  snow  avalanches  one  must  remain  cognizant  of 
the  fact  that  the  cohesive  nature  of  snow  has  been  neglected.  This  omission  is  probably  most 
'mpractical  for  modeling  the  initiation  of  the  avalanche  and  the  final  stages  of  compaction  in 
the  runout  zone,  but  it  is  unknown  if  neglecting  cohesion  plays  a  significant  role  in  determining 
the  salient  features  of  the  flow.  The  same  statement  may  be  made  concerning  particle  size 
distribution.  It  is  sufficient  to  say  that  some  progress  has  been  made  in  modeling  granular 
flows,  but  efforts  should  continue  to  resolve  these  and  many  other  issues. 
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